﻿using MathNet.Numerics;
using MathNet.Numerics.LinearAlgebra;
using System;
using System.Collections;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.Windows.Forms;

namespace CurveFitting
{
    public partial class Form1 : Form
    {
        public Form1()
        {
            InitializeComponent();
            InitPanel();
        }
        public void InitPanel()
        {
            double[] xValues = { 0.251, 0.28, 0.364, 0.573, 0.588, 0.599 };
            double[] yValues = { 0, 2.5, 5, 10, 20, 40 };

            int dimensionl = 3;
            double[] result = CurveFitHelper.MultiLine(xValues, yValues, xValues.Length, dimensionl);

            //计算出拟合曲线后的第一个Y值
            double yFristValue = 0;
            for (int j = 0; j < dimensionl + 1; j++)
            {
                yFristValue += result[j] * Math.Pow(xValues[0], j);
            }

            ArrayList stationaryPoint = CurveFitHelper.GetStationaryPoint(result, 100, xValues[0], xValues[xValues.Length - 1], yFristValue);

            stationaryPoint.Insert(0, xValues[0]);
            stationaryPoint.Insert(stationaryPoint.Count, xValues[xValues.Length - 1]);

            float[] yy = new float[stationaryPoint.Count];
            string[] xx = new string[stationaryPoint.Count];
            int count = 0;
            foreach (double values in stationaryPoint)
            {
                xx[count] = values.ToString("F3");
                double sum = 0;
                for (int j = 0; j < dimensionl + 1; j++)
                {
                    sum += result[j] * Math.Pow(values, j);
                }
                yy[count] = (float)sum;
                count++;
            }
            panel1.BackgroundImage = Curve2D.DrawImage(panel1.Width, panel1.Height, xx, yy);

        }
        private void button1_Click(object sender, EventArgs e)
        {
            double[] xValues = { 0.251, 0.28, 0.364, 0.573, 0.588, 0.599 };
            double[] yValues = { 0, 2.5, 5, 10, 20, 40 };
            Random rand = new Random();
            //for (int i = 0; i < 10; i++)
            //{
            //    xValues[i] = rand.Next(20, 40);
            //    yValues[i] = rand.Next(20, 40);
            //}
            //Array.Sort(xValues);
            //Array.Sort(yValues);
            //Array.Reverse(xValues);
            //Array.Reverse(yValues);

            double[] res = FittingCurveByLeastSquare(xValues, yValues);
            double[] res1 = CurveFitHelper.MultiLineToFixedPoint(xValues, yValues, xValues.Length, 2, 0, 0);
            double[] res2 = CurveFitHelper.MultiLineToFixedPoint2(xValues, yValues, xValues.Length, 2, 0, 0);
            double[] res3 = CurveFitHelper.MultiLineToOrigin(xValues, yValues, xValues.Length, 2);

            double[] yy = new double[yValues.Length];
            for (int i = 0; i < xValues.Length; i++)
            {
                double result = 0;
                for (int j = 0; j < res2.Length; j++)
                {
                    result += res2[j] * Math.Pow(xValues[i], j);
                }
                yy[i] = result;
            }
            double R = CurveFitHelper.CalculateRSquared(yValues, yy);
            this.textBox1.Text = res[0].ToString("F2") + "x^2" + " + " + res[1].ToString("F2") + "x" + " + " + res[2].ToString("F2");
        }
        /// <summary>
        /// 拟合曲线所筛选的点的个数
        /// </summary>
        const int CurveNums = 10;

        /// <summary>
        /// 拟合曲线的阶数
        /// </summary>
        public int m = 2;

        /// <summary>
        /// 取直线上的点的个数
        /// </summary>
        const int XNums = 50;
        /// <summary>
        ///  最小二乘法拟合曲线
        /// </summary>
        /// <param name="X">X轴数组</param>
        /// <param name="Y">Y轴数组</param>
        /// <param name="m">阶数</param>
        /// <returns>返回曲线方程的各阶系数（由高阶到低阶,一般m=2)</returns>
        public double[] FittingCurveByLeastSquare(double[] X, double[] Y, int m = 2)
        {
            ///
            double[] res = new double[m + 1];
            if (X.Length > m && Y.Length > m)
            {
                res = Fit.Polynomial(X, Y, m);
            }
            return res;
        }

        /// <summary>
        /// 计算直线方程
        /// </summary>
        /// <param name="StartPoint">直线起点</param>
        /// <param name="angle">直线的角度</param>
        /// <returns>返回k,b的数组</returns>
        public double[] CalculateLine(Point StartPoint, double angle)
        {
            double k = Math.Tan(angle / 180 * Math.PI);
            double b = StartPoint.Y - k * StartPoint.X;
            return new double[] { k, b };
        }

        /// <summary>
        /// 计算交点
        /// </summary>
        /// <param name="LineResult">直线的k,b</param>
        /// <param name="CurveResult">拟合曲线的系数数组</param>
        /// <param name="LineX">直线上的点的X坐标</param>
        /// <param name="m">阶数</param>
        /// <returns></returns>
        public Point CalculateInterPoint(double[] LineResult, double[] CurveResult, double[] LineX, int m)
        {
            Point InterPoint = new Point();
            if (LineResult.Length == 2 && CurveResult.Length == (m + 1))
            {
                double k = LineResult[0];
                double b = LineResult[1];
                for (int i = 0; i < LineX.Count(); i++)
                {
                    double x = LineX[i];
                    double y_Line = k * x + b;
                    double y_Curve = 0;
                    for (int n = 0; n <= m; n++)
                    {
                        y_Curve += CurveResult[n] * Math.Pow(x, n);
                    }
                    bool IsSuccessFind = false;
                    for (int t = 1; t <= 50; t++)
                    {
                        if (Math.Abs(y_Line - y_Curve) < 5)
                        {
                            InterPoint = new Point((int)x, (int)y_Line);
                            IsSuccessFind = true;
                            break;//跳出内循环
                        }
                    }
                    if (IsSuccessFind)
                    {
                        break;//跳出外循环
                    }
                }
            }
            return InterPoint;
        }

        /// <summary>
        /// 计算R^2,R^2这个值越接近1，说明拟合出来的曲线跟原曲线就越接近
        /// </summary>
        /// <param name="Y">实际的Y</param>
        /// <param name="Ytest">代入拟合曲线方程得到的Y</param>
        /// <returns>返回R^2</returns>
        public double CalculateRSquared(double[] Y, double[] Ytest)
        {
            double RSquared = GoodnessOfFit.RSquared(Y, Ytest);
            return RSquared;
        }
    }
}
